﻿

namespace LxCFD
{
    public class ConnectWallOneSide : ConnectWallBase
    {
        public void setConnect(ValueBase<double> S, EquationBase OMatrix, EquationBase NMatrix)
        {
            setArray_2(S, OMatrix, NMatrix);
        }
        public override void SetAb_Step1(int taskIndex)
        {
            if (S is Value1Base<double> s1)
            {
                if (s1.IsEqual0()) 
                {
                    InnerIndex.FaceIDs.GetSeg(out int i1, out int endp1, taskIndex);
                    for (; i1 < endp1; i1++)
                    {
                        ref PairMatrix ab = ref ADLUb[i1];
                        ab.b2 = 0;
                        ab.a22 = 0;

                        int ofid = i1 + deltaOFaceId;
                        int nfid = i1 + deltaNFaceId;
                        setAb_oneSideWall(ref ab, ofid, ref O_Faces[ofid], nfid, ref N_Faces[nfid]);
                    }
                }
                else
                {
                    InnerIndex.FaceIDs.GetSeg(out int i1, out int endp1, taskIndex);
                    for (; i1 < endp1; i1++)
                    {
                        ref PairMatrix ab = ref ADLUb[i1];
                        int ofid = i1 + deltaOFaceId;
                        ref FvFace oface = ref O_Faces[ofid];
                        ab.b2 = s1.GetFaceValue(ofid) * oface.Area;
                        ab.a22 = 0;

                        int nfid = i1 + deltaNFaceId;
                        setAb_oneSideWall(ref ab, ofid, ref oface, nfid, ref N_Faces[nfid]);
                    }
                }
            }
            else if (S is Value2Base<double> s2)
            {
                InnerIndex.FaceIDs.GetSeg(out int i1, out int endp1, taskIndex);
                for (; i1 < endp1; i1++)
                {
                    ref PairMatrix ab = ref ADLUb[i1];
                    int ofid = i1 + deltaOFaceId;
                    ref FvFace oface = ref O_Faces[ofid];
                    s2.GetFaceValue(out double sou_ad, out double sou_b, ofid);
                    ab.b2 = sou_b * oface.Area;
                    ab.a22 = -sou_ad * oface.Area;

                    int nfid = i1 + deltaNFaceId;
                    setAb_oneSideWall(ref ab, ofid, ref oface, nfid, ref N_Faces[nfid]);
                }
            }
            else
            {
                LxConsole.Error(string.Format($"接触面{Contact.name}源项无法识别"));
            }
        }
        protected void setAb_oneSideWall(ref PairMatrix ab, int ofid, ref FvFace oface, int nfid, ref FvFace nface)
        {
            double ogamma = O_Gamma.GetFaceValue(ofid);
            double ogg = ogamma * oface.gDiff;
            double ofacecor = ogamma * oface.CalNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_Cells); 
            ocor_a *= ogamma;
            ocor_b *= ogamma;

            double ngamma = N_Gamma.GetFaceValue(nfid);
            double ngg = ngamma * nface.gDiff;
            double nfacecor = ngamma * nface.CalNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_Cells);
            ncor_a *= ngamma;
            ncor_b *= ngamma;

            ab.a11 = ogg; ab.a12 = -ogg; ab.b1 = ofacecor;
            ab.a44 = ngg; ab.a43 = -ngg; ab.b4 = nfacecor;
            ab.a22 += ogg + ngg + ocor_a + ncor_a;//一次成型
            ab.a21 = -ogg - ocor_a;
            ab.a24 = -ngg - ncor_a;
            ab.b2 += -(ocor_b + ncor_b);
            //第三行均没有赋值，也没有初始化,直接按等于第二行赋值
            ab.a31 = ab.a21;
            ab.a33 = ab.a22;
            ab.a34 = ab.a24;
            ab.b3 = ab.b2;
            ab.Eliminate_a12a43();
            setAb_ONALU(ref ab, ofid, nfid);
        }
        public override void Update_x_Bound(int taskIndex)
        {
            InnerIndex.FaceIDs.GetSeg(out int i1, out int endp1, taskIndex);
            for (; i1 < endp1; i1++)
            {
                ref FvPair pair = ref Pairs[i1];
                ref PairMatrix ab = ref ADLUb[i1];
                double x1 = O_x[pair.O_InnCell];
                double x4 = N_x[pair.N_InnCell];
                O_x[pair.O_Cell] = N_x[pair.N_Cell] = (ab.b2 - ab.a21 * x1 - ab.a24 * x4) / ab.a22;
            }
        }
    }
}
//iniSource(taskIndex);
//InnerIndex.FaceIDs.GetSeg(out int i1, out int endp1, taskIndex);
//for (; i1 < endp1; i1++)
//{
//    int ofid = i1 + deltaOFaceId;
//    ref Face oface = ref O_faces[ofid];
//    double ogamma = O_Gamma.getFaceValue(ofid);
//    double ogg = ogamma * oface.gDiff;
//    double ofacecor = ogamma * oface.CalNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_cells);
//    ocor_a *= ogamma;
//    ocor_b *= ogamma;

//    int nfid = i1 + deltaNFaceId;
//    ref Face nface = ref N_faces[nfid];
//    double ngamma = N_Gamma.getFaceValue(nfid);
//    double ngg = ngamma * nface.gDiff;
//    double nfacecor = ngamma * nface.CalNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_cells);
//    ncor_a *= ngamma;
//    ncor_b *= ngamma;

//    ref PairMatrixAb ab = ref ADLUb[i1];
//    ab.a11 = ogg; ab.a12 = -ogg; ab.b1 = ofacecor;
//    ab.a44 = ngg; ab.a43 = -ngg; ab.b4 = nfacecor;
//    ab.a22 += ogg + ngg + ocor_a + ncor_a;//一次成型
//    ab.a21 = -ogg - ocor_a;
//    ab.a24 = -ngg - ncor_a;
//    ab.b2 += -(ocor_b + ncor_b);
//    //第三行均没有赋值
//    ab.a31 = ab.a21;
//    ab.a33 = ab.a22;
//    ab.a34 = ab.a24;
//    ab.b3 = ab.b2;
//    ab.Eliminate_a12a43();
//    O_ALUN[ofid] = -ab.a11;
//    O_ALUO[ofid] = 1;
//    O_bFaceCommon[ofid] = ab.b1;
//    N_ALUN[nfid] = -ab.a44;
//    N_ALUO[nfid] = 1;
//    N_bFaceCommon[nfid] = ab.b4;
//}
#region 备份
//InnerIndex.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//for (; i1 < endp1; i1++)
//{
//    int ofid = i1 + deltaOFaceId;
//    ref Face oface = ref O_faces[ofid];
//    double ogamma = O_Gamma.getFaceValue(ofid);
//    double ogg = ogamma * oface.gDiff;
//    double ofacecor = ogamma * oface.calNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_cells);
//    ocor_a *= ogamma;
//    ocor_b *= ogamma;

//    int nfid = i1 + deltaNFaceId;
//    ref Face nface = ref N_faces[nfid];
//    double ngamma = N_Gamma.getFaceValue(nfid);
//    double ngg = ngamma * nface.gDiff;
//    double nfacecor = ngamma * nface.calNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_cells);
//    ncor_a *= ngamma;
//    ncor_b *= ngamma;

//    ref PairMatrixAb ab = ref ADLUb[i1];
//    ab.a12 = ab.a21 = -ogg;
//    ab.a43 = ab.a24 = -ngg;
//    ab.a11 = ogg; ab.b1 = ofacecor;
//    ab.a44 = ngg; ab.b4 = nfacecor;
//    ab.a22 += ogg + ngg + ocor_a + ncor_a;//一次成型
//    ab.b2 += -(ocor_b + ncor_b);
//    //第三行均没有赋值
//    ab.a31 = ab.a21;
//    ab.a33 = ab.a22;
//    ab.a34 = ab.a24;
//    ab.b3 = ab.b2;
//    ab.Eliminate_a12a43();
//    O_ALUN[ofid] = -ab.a11;
//    O_ALUO[ofid] = 1;
//    O_bFaceCommon[ofid] = ab.b1;
//    N_ALUN[nfid] = -ab.a44;
//    N_ALUO[nfid] = 1;
//    N_bFaceCommon[nfid] = ab.b4;
//}
#endregion
#region 暂时封存
//InnerIndex.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//for (; i1 < endp1; i1++)
//{
//    int ofid = i1 + deltaOFaceId;
//    ref Face oface = ref O_faces[ofid];
//    double ogamma = O_Gamma.getFaceValue(ofid);
//    double ogg = ogamma * oface.gDiff;
//    double ofacecor = ogamma * oface.calNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_cells);
//    ocor_a *= ogamma;
//    ocor_b *= ogamma;
//    O_ALUO[ofid] = O_ALUN[ofid] = -ogg;
//    O_bFaceCommon[ofid] = ofacecor;

//    int nfid = i1 + deltaNFaceId;
//    ref Face nface = ref N_faces[nfid];
//    double ngamma = N_Gamma.getFaceValue(nfid);
//    double ngg = ngamma * nface.gDiff;
//    double nfacecor = ngamma * nface.calNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_cells);
//    ncor_a *= ngamma;
//    ncor_b *= ngamma;
//    N_ALUO[nfid] = N_ALUN[nfid] = -ngg;
//    N_bFaceCommon[nfid] = nfacecor;

//    O_ADiag[oface.NeighborCellIndex] = ogg + ngg + ocor_a + ncor_a;//一次成型
//    O_b[oface.NeighborCellIndex] = -(ocor_b + ncor_b);
//}
#endregion
#region BiCGStab
//public override void Initialize_Step1_r0rp1rho(int taskIndex)
//{
//    //var r0 = BiCG.O_Solver.r0;
//    //var r = BiCG.O_Solver.r;
//    //var p = BiCG.O_Solver.p;
//    //InnerIndex.CellIdx.GetSeg(out int i1, out int endp1, taskIndex);
//    //double temp = 0;
//    //for (; i1 < endp1; i1++)
//    //{
//    //    ref Mesh.Pair pair = ref Pairs[i1];
//    //    double rt = r[pair.O_Cell] = r0[pair.O_Cell] = p[pair.O_Cell]
//    //        = O_b[pair.O_Cell] - calALUProductValue(ref pair, O_x, N_x) - O_ADiag[pair.O_Cell] * O_x[pair.O_Cell];
//    //    temp += rt * rt;
//    //}
//    //TasksControl.Storage1[taskIndex] += temp;
//}
//public override void Cal_rhoPart1(int taskIndex)
//{
//    BiCGBase.DotProduct(O_Index, BiCG.O_Solver.r0, BiCG.O_Solver.r, taskIndex);
//}
//public override void Cal_pArray(int taskIndex)
//{
//    BiCGBase.Cal_pArray(O_Index, BiCG.O_Solver.p, BiCG.O_Solver.r, BiCG.O_Solver.v, taskIndex);
//}
//public override void Cal_yArray(int taskIndex)
//{
//    BiCGBase.DInversePrecond(O_Index, BiCG.O_Solver.y, O_ADiag, BiCG.O_Solver.p, taskIndex);
//}
//public override void Copy_yArray(int taskIndex)
//{
//    copy_Array(BiCG.O_Solver.y, BiCG.N_Solver.y, taskIndex);
//}
//public override void Cal_vArray_alpha1(int taskIndex)
//{
//    double[] v = BiCG.O_Solver.v;
//    double[] y = BiCG.O_Solver.y;
//    double[] Ny = BiCG.N_Solver.y;
//    double[] r0 = BiCG.O_Solver.r0;
//    InnerIndex.CellIdx.GetSeg(out int i1, out int endp1, taskIndex);
//    double temp = 0;
//    for (; i1 < endp1; i1++)
//    {
//        ref Mesh.Pair pair = ref Pairs[i1];
//        double vt = v[pair.O_Cell] = O_ADiag[pair.O_Cell] * y[pair.O_Cell] + calALUProductValue(ref pair, y, Ny);
//        temp += r0[pair.O_Cell] * vt;
//    }
//    TasksControl.Storage1[taskIndex] += temp;
//}
//public override void Cal_sArray(int taskIndex)
//{
//    BiCGBase.Cal_rMinus(O_Index, BiCG.O_Solver.r, BiCGBase.alpha, BiCG.O_Solver.v, taskIndex);
//}
//public override void Cal_zArray(int taskIndex)
//{
//    BiCGBase.DInversePrecond(O_Index, BiCG.O_Solver.z, O_ADiag, BiCG.O_Solver.r, taskIndex);
//}
//public override void Copy_zArray(int taskIndex)
//{
//    copy_Array(BiCG.O_Solver.z, BiCG.N_Solver.z, taskIndex);
//}
//public override void Cal_tArray_omega1(int taskIndex)
//{
//    double[] s = BiCG.O_Solver.r;
//    double[] z = BiCG.O_Solver.z;
//    double[] Nz = BiCG.N_Solver.z;
//    double[] t = BiCG.O_Solver.t;
//    InnerIndex.CellIdx.GetSeg(out int i1, out int endp1, taskIndex);
//    double temp1 = 0, temp2 = 0;
//    for (; i1 < endp1; i1++)
//    {
//        ref Mesh.Pair pair = ref Pairs[i1];
//        double st = s[pair.O_Cell];
//        double tt = t[pair.O_Cell]
//            = O_ADiag[pair.O_Cell] * z[pair.O_Cell]
//            + calALUProductValue(ref pair, z, Nz);
//        temp1 += tt * st;
//        temp2 += tt * tt;
//    }
//    TasksControl.Storage1[taskIndex] += temp1;
//    TasksControl.Storage2[taskIndex] += temp2;
//}
//public override void Cal_xArray(int taskIndex)
//{
//    BiCGBase.cal_xArray(O_Index, O_x, BiCG.O_Solver.y, BiCG.O_Solver.z, taskIndex);
//}
//protected void copy_Array(double[] FromArray, double[] ToArray, int taskIndex)
//{
//    InnerIndex.CellIdx.GetSeg(out int i1, out int endp1, taskIndex);
//    for (; i1 < endp1; i1++)
//    {
//        ref Mesh.Pair pair = ref Pairs[i1];
//        ToArray[pair.N_Cell] = FromArray[pair.O_Cell];
//    }
//}
#endregion
//public override void Relax(int taskIndex)
//{
//    if (relaxFactor == 1) return;
//    double relax1 = (1 - relaxFactor);
//    O_Index.CellIdx.GetSeg(out int i1, out int endp1, taskIndex);
//    for (; i1 < endp1; i1++)
//    {
//        ref double adiag = ref O_ADiag[i1];
//        adiag /= relaxFactor;
//        O_b[i1] += relax1 * adiag * O_x[i1];
//    }
//}
//public override void CalByGS(int taskIndex)
//{
//    cal_x(O_x, N_x, O_x, N_x, taskIndex);
//}
//public override void CalByJacobi_Step1(int taskIndex)
//{
//    cal_x(O_x, N_x, O_xCommon, N_xCommon, taskIndex);
//}
//public override void CalByJacobi_Step2(int taskIndex)
//{
//    cal_x(O_xCommon, N_xCommon, O_x, N_x, taskIndex);
//}
//protected void cal_x(double[] O_from, double[] N_from, double[] O_to, double[] N_to, int taskIndex)
//{
//    InnerIndex.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//    for (; i1 < endp1; i1++)
//    {
//        ref Mesh.Pair pair = ref Pairs[i1];
//        N_to[pair.N_Cell] = O_to[pair.O_Cell] = (O_b[pair.O_Cell]
//            - O_ALUN[pair.O_Face] * O_from[pair.O_InnCell]
//            - N_ALUN[pair.N_Face] * N_from[pair.N_InnCell]
//              ) / O_ADiag[pair.O_Cell];
//    }
//}
//protected double calALUProductValue(ref Mesh.Pair pair, double[] OFromArray, double[] NFromArray)
//{
//    return O_ALUN[pair.O_Face] * OFromArray[pair.O_InnCell]
//            + N_ALUN[pair.N_Face] * NFromArray[pair.N_InnCell];
//}
//protected double cal_rMinusALU(ref Mesh.Pair pair)
//{
//    return O_b[pair.O_Cell]
//            - O_ALUN[pair.O_Face] * O_x[pair.O_InnCell]
//            - N_ALUN[pair.N_Face] * N_x[pair.N_InnCell];
//}
#region 接口IPassDataToCoarserContact
//public override void SendRestrict_ALUr(out double A_Diag, out double OZone_A_LU_O, out double OZone_A_LU_N, out double NZone_A_LU_O, out double NZone_A_LU_N, ref MGMesh.Pair CoarsePair)
//{
//    double atemp = 0, rtemp = 0, aluo = 0, alun = 0;
//    foreach (int fpair in CoarsePair.FinePairs)
//    {
//        ref Mesh.Pair pair = ref Pairs[fpair];
//        double atm = O_ADiag[pair.O_Cell];
//        atemp += atm;

//        double ozalu = O_ALUN[pair.O_Face];
//        aluo += ozalu;

//        double nzalu = N_ALUN[pair.N_Face];
//        alun += nzalu;

//        rtemp += O_b[pair.O_Cell]
//            - ozalu * O_x[pair.O_InnCell]
//            - nzalu * N_x[pair.N_InnCell]
//            - atm * O_x[pair.O_Cell];
//    }
//    A_Diag = atemp;
//    CoarsePair.r = rtemp;
//    OZone_A_LU_O = OZone_A_LU_N = aluo;
//    NZone_A_LU_O = NZone_A_LU_N = alun;
//}
//public override void SendRestrict_r(ref MGMesh.Pair CoarsePair)
//{
//    double rtemp = 0;
//    foreach (int fpair in CoarsePair.FinePairs)
//    {
//        rtemp += cal_rMinusALU(ref Pairs[fpair]);
//    }
//    CoarsePair.r = rtemp;
//}
#endregion
//public void Pass_A_r_ToCoarser(out double A_owner, ref MGMesh.Cell ownerCoarser)
//{
//    double Atemp = 0;
//    double rtemp = 0;
//    foreach (int ocell in ownerCoarser.FineCells)
//    {
//        int oface = O_cells[ocell].CellFacesAsN[0];
//        int nface = oface + deltaFaceIndex;
//        int ncell = ocell + deltaCellIndex;

//        double atemp = O_A_Diag[ocell];
//        Atemp += atemp;
//        rtemp += O_b[ocell]
//            - O_A_LU[oface] * O_x[O_faces[oface].OwnerCellIndex]
//            - N_A_LU[nface] * N_x[N_faces[nface].OwnerCellIndex]
//            - atemp * O_x[ocell];
//    }
//    A_owner = Atemp;
//    ownerCoarser.r = rtemp;
//}
//public void Pass_r_ToCoarser(ref MGMesh.Cell ownerCoarser)
//{
//    double rtemp = 0;
//    foreach (int ocell in ownerCoarser.FineCells)
//    {
//        int oface = O_cells[ocell].CellFacesAsN[0];
//        int nface = oface + deltaFaceIndex;
//        int ncell = ocell + deltaCellIndex;

//        double atemp = O_A_Diag[ocell];
//        rtemp += O_b[ocell]
//            - O_A_LU[oface] * O_x[O_faces[oface].OwnerCellIndex]
//            - N_A_LU[nface] * N_x[N_faces[nface].OwnerCellIndex]
//            - atemp * O_x[ocell];
//    }
//    ownerCoarser.r = rtemp;
//}
//public class ConnectOfOneSideWall : ConnectBase, IConnect
//{
//    protected double relaxFactor = 1;
//    public void setConnect(Value12Base S)
//    {
//        setArray(S);

//        this.relaxFactor = Math.Min(ownerScalar.relaxFactor, neighborScalar.relaxFactor);
//    }

//    public void SetAb_Step1(int taskIndex)
//    {
//        InnerIndex.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//        for (; i1 < endp1; i1++)
//        {
//            int ofid = i1 + deltaOFaceId;
//            ref Face oface = ref O_faces[ofid];
//            double ogamma = O_Gamma.getFaceValue(ofid);
//            double ogg = ogamma * oface.gDiff;
//            double ofacecor = ogamma * oface.calNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_cells);
//            ocor_a *= ogamma;
//            ocor_b *= ogamma;
//            O_A_LU_N[ofid] = -ogg;
//            O_bFaceCommon[ofid] = ofacecor;

//            int nfid = i1 + deltaNFaceId;
//            ref Face nface = ref N_faces[nfid];
//            double ngamma = N_Gamma.getFaceValue(nfid);
//            double ngg = ngamma * nface.gDiff;
//            double nfacecor = ngamma * nface.calNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_cells);
//            ncor_a *= ngamma;
//            ncor_b *= ngamma;
//            N_A_LU_N[nfid] = -ngg;
//            N_bFaceCommon[nfid] = nfacecor;

//            O_A_Diag[oface.NeighborCellIndex] = ogg + ngg + ocor_a + ncor_a;//一次成型
//            O_b[oface.NeighborCellIndex] = -(ocor_b + ncor_b);
//        }
//        ///以下代码暂时不要删，否则反应不过来了
//        //O_Index.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//        //for (; i1 < endp1; i1++)
//        //{
//        //    ref Face oface = ref O_faces[i1];
//        //    double ogamma = O_Gamma.getFaceValue(i1);
//        //    double ogg = ogamma * oface.gDiff;
//        //    double ofacecor = ogamma * oface.calNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_cells);
//        //    ocor_a *= ogamma;
//        //    ocor_b *= ogamma;
//        //    O_A_LU_N[i1] = -ogg;
//        //    O_bFaceCommon[i1] = ofacecor;

//        //    int nfidx = i1 + deltaFaceIndex;
//        //    ref Face nface = ref N_faces[nfidx];
//        //    double ngamma = N_Gamma.getFaceValue(nfidx);
//        //    double ngg = ngamma * nface.gDiff;
//        //    double nfacecor = ngamma * nface.calNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_cells);
//        //    ncor_a *= ngamma;
//        //    ncor_b *= ngamma;
//        //    N_A_LU_N[nfidx] = -ngg;
//        //    N_bFaceCommon[nfidx] = nfacecor;

//        //    O_A_Diag[oface.NeighborCellIndex] = ogg + ngg + ocor_a + ncor_a;//一次成型
//        //    O_b[oface.NeighborCellIndex] = -(ocor_b + ncor_b);
//        //}
//    }
//    public void SetAb_Step2(int taskIndex)
//    {
//        if (S is Value1Base s1)
//        {
//            if (s1 is Value1Const s1c && s1c.number == 0) return;
//            O_Index.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//            for (; i1 < endp1; i1++)
//            {
//                ref Face face = ref O_faces[i1];
//                O_b[face.NeighborCellIndex] += s1.getFaceValue(i1) * face.Area;
//            }

//        }
//        else if (S is Value2Base s2)
//        {
//            O_Index.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//            for (; i1 < endp1; i1++)
//            {
//                ref Face face = ref O_faces[i1];
//                s2.getFaceValue(out double sou_ad, out double sou_b, i1);
//                O_A_Diag[face.NeighborCellIndex] -= sou_ad * face.Area;
//                O_b[face.NeighborCellIndex] += sou_b * face.Area;
//            }
//        }
//        else
//        {
//            LxConsole.AddLine(LogTypes.错误, string.Format($"接触面{contact.name}源项无法识别"));
//        }
//    }


//    public void Relax(int taskIndex)
//    {
//        double relax1 = (1 - relaxFactor);
//        O_Index.CellIdx.GetSeg(out int i1, out int endp1, taskIndex);
//        for (; i1 < endp1; i1++)
//        {
//            O_A_Diag[i1] /= relaxFactor;
//            O_b[i1] += relax1 * O_A_Diag[i1] * O_x[i1];
//        }
//    }
//    //public void Cal_x_GS(int taskIndex)
//    //{
//    //    cal_x(O_x, N_x, O_x, N_x, taskIndex);
//    //}
//    //public void Cal_x_JacobiStep1(int taskIndex)
//    //{
//    //    cal_x(O_x, N_x, O_xCommon, N_xCommon, taskIndex);
//    //}
//    //public void Cal_x_JacobiStep2(int taskIndex)
//    //{
//    //    cal_x(O_xCommon, N_xCommon, O_x, N_x, taskIndex);
//    //}
//    //protected void cal_x(double[] O_from, double[] N_from, double[] O_to, double[] N_to, int taskIndex)
//    //{
//    //    O_Index.FaceIdx.GetSeg(out int oface, out int endp1, taskIndex);
//    //    for (; oface < endp1; oface++)
//    //    {
//    //        int ocell = O_faces[oface].NeighborCellIndex;
//    //        int nface = oface + NFaceBaseO;
//    //        int ncell = ocell + NCellBaseO;

//    //        N_to[ncell] = O_to[ocell] = (O_b[ocell]
//    //            - O_A_LU_N[oface] * O_from[O_faces[oface].OwnerCellIndex]
//    //            - N_A_LU_N[nface] * N_from[N_faces[nface].OwnerCellIndex]
//    //              ) / O_A_Diag[ocell];
//    //    }
//    //}
//    #region 接口IPassDataToCoarserContact
//    //public void pass_A_r_owner_ToCoarser(out double A_owner, ref MGMesh.Cell ownerCoarser)
//    //{
//    //    double Atemp = 0;
//    //    double rtemp = 0;
//    //    foreach (int ocell in ownerCoarser.FineCells)
//    //    {
//    //        int oface = O_cells[ocell].CellFacesAsN[0];
//    //        int nface = oface + NFaceBaseO;
//    //        int ncell = ocell + NCellBaseO;

//    //        double atemp = O_A_Diag[ocell];
//    //        Atemp += atemp;
//    //        rtemp += O_b[ocell]
//    //            - O_A_LU_N[oface] * O_x[O_faces[oface].OwnerCellIndex]
//    //            - N_A_LU_N[nface] * N_x[N_faces[nface].OwnerCellIndex]
//    //            - atemp * O_x[ocell];
//    //    }
//    //    A_owner = Atemp;
//    //    ownerCoarser.r = rtemp;
//    //}
//    //public void pass_r_owner_ToCoarser(ref MGMesh.Cell ownerCoarser)
//    //{
//    //    double rtemp = 0;
//    //    foreach (int ocell in ownerCoarser.FineCells)
//    //    {
//    //        int oface = O_cells[ocell].CellFacesAsN[0];
//    //        int nface = oface + NFaceBaseO;
//    //        int ncell = ocell + NCellBaseO;

//    //        double atemp = O_A_Diag[ocell];
//    //        rtemp += O_b[ocell]
//    //            - O_A_LU_N[oface] * O_x[O_faces[oface].OwnerCellIndex]
//    //            - N_A_LU_N[nface] * N_x[N_faces[nface].OwnerCellIndex]
//    //            - atemp * O_x[ocell];
//    //    }
//    //    ownerCoarser.r = rtemp;
//    //}

//    #endregion
//    public struct MatrixToEliminateBound
//    {
//        public A3b1 r1;
//        public A3b1 r2;
//        public A3b1 r3;
//        public void Eliminate()
//        {
//            r1.Minus(r1.a2 / r2.a2, r2);
//            r3.Minus(r3.a2 / r2.a2, r2);
//        }
//        public void Print()
//        {
//            Console.WriteLine(r1);
//            Console.WriteLine(r2);
//            Console.WriteLine(r3);
//        }
//    }
//    public struct A3b1
//    {
//        public double a1;
//        public double a2;
//        public double a3;
//        public double b;
//        public void Set(double a1, double a2, double a3, double b)
//        {
//            this.a1 = a1;
//            this.a2 = a2;
//            this.a3 = a3;
//            this.b = b;
//        }
//        public void Multiply(double d)
//        {
//            a1 *= d;
//            a2 *= d;
//            a3 *= d;
//            b *= d;
//        }
//        public void Minus(double d, A3b1 ab)
//        {
//            a1 -= d * ab.a1;
//            a2 -= d * ab.a2;
//            a3 -= d * ab.a3;
//            b -= d * ab.b;
//        }
//        public override string ToString()
//        {
//            return string.Format($"| {a1} | {a2} | {a3} | {b} |");
//        }
//    }
//}

//protected void SetAb_Step2(int taskIndex)
//{
//    if (S is Value1Base s1)
//    {
//        if (s1 is Value1Const s1c && s1c.number == 0)
//        {
//            return;
//            //InnerIndex.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//            //for (; i1 < endp1; i1++)
//            //{
//            //    //ref Face face = ref O_faces[i1];
//            //    //O_A_Diag[face.NeighborCellIndex] = 0;
//            //    //O_b[face.NeighborCellIndex] = 0;
//            //    Local_A_Diag[i1] = 0;
//            //    Local_b[i1] = 0;
//            //}
//        }
//        else
//        {
//            InnerIndex.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//            for (; i1 < endp1; i1++)
//            {
//                int ofaceid = i1 + deltaOFaceId;
//                ref Face face = ref O_faces[ofaceid];
//                Local_A_Diag[i1] = 0;
//                Local_b[i1] = s1.getFaceValue(ofaceid) * face.Area;
//                //O_A_Diag[face.NeighborCellIndex] = 0;
//                //O_b[face.NeighborCellIndex] = s1.getFaceValue(i1) * face.Area;
//            }
//        }
//    }
//    else if (S is Value2Base s2)
//    {
//        InnerIndex.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//        for (; i1 < endp1; i1++)
//        {
//            //ref Face face = ref O_faces[i1];
//            //s2.getFaceValue(out double sou_ad, out double sou_b, i1);
//            //O_A_Diag[face.NeighborCellIndex] = -sou_ad * face.Area;
//            //O_b[face.NeighborCellIndex] = sou_b * face.Area;
//            int ofaceid = i1 + deltaOFaceId;
//            ref Face face = ref O_faces[ofaceid];
//            s2.getFaceValue(out double sou_ad, out double sou_b, ofaceid);
//            Local_A_Diag[i1] = -sou_ad * face.Area;
//            Local_b[i1] = sou_b * face.Area;
//        }
//    }
//    else
//    {
//        LxConsole.AddLine(LogTypes.错误, string.Format($"接触面{contact.name}源项无法识别"));
//    }
//}
//protected void setAb(int iid, int ofid, ref Face oface, double AbySource, double bbySource)
//{
//    int nfid = iid + deltaNFaceId;
//    ref double alun = ref Local_A_LU_N[iid];
//    ref double aluo = ref Local_A_LU_O[iid];

//    //ref Face oface = ref O_faces[ofid];
//    double ogamma = O_Gamma.getFaceValue(ofid);
//    double ogg = ogamma * oface.gDiff;
//    double ofacecor = ogamma * oface.calNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_cells);
//    ocor_a *= ogamma;
//    ocor_b *= ogamma;
//    aluo = -ogg;

//    ref Face nface = ref N_faces[nfid];
//    double ngamma = N_Gamma.getFaceValue(nfid);
//    double ngg = ngamma * nface.gDiff;
//    double nfacecor = ngamma * nface.calNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_cells);
//    ncor_a *= ngamma;
//    ncor_b *= ngamma;
//    alun = -ngg;

//    ref double adiag = ref Local_A_Diag[iid];
//    ref double b = ref Local_b[iid];
//    adiag = AbySource + ogg + ngg + ocor_a + ncor_a;
//    b = bbySource - (ocor_b + ncor_b);

//    MatrixToEliminateBound mtool = new MatrixToEliminateBound();
//    mtool.r1.Set(-aluo, aluo, 0, ofacecor);
//    mtool.r2.Set(aluo, adiag, alun, b);
//    mtool.r3.Set(0, alun, -alun, nfacecor);
//    mtool.Eliminate();
//    O_A_LU_N[ofid] = -mtool.r1.a1;//用以计算A_Diag，还得取反
//    N_A_LU_N[nfid] = -mtool.r3.a3;
//    BetZone_A_LU_O[iid] = mtool.r1.a3;
//    BetZone_A_LU_N[iid] = mtool.r3.a1;
//    O_bFaceCommon[ofid] = mtool.r1.b;
//    N_bFaceCommon[nfid] = mtool.r3.b;

//    O_A_LU_O[ofid] = N_A_LU_O[nfid] = 1;//这是与对面区域的系数
//}
//public void SetAb_Step1(int taskIndex)
//{
//    InnerIndex.FaceIdx.GetSeg(out int iid, out int endp1, taskIndex);
//    MatrixToEliminateBound mtool = new MatrixToEliminateBound();
//    for (; iid < endp1; iid++)
//    {
//        int ofid = iid + deltaOFaceId;
//        int nfid = iid + deltaNFaceId;
//        ref double alun = ref Local_A_LU_N[iid];
//        ref double aluo = ref Local_A_LU_O[iid];

//        ref Face oface = ref O_faces[ofid];
//        double ogamma = O_Gamma.getFaceValue(ofid);
//        double ogg = ogamma * oface.gDiff;
//        double ofacecor = ogamma * oface.calNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_cells);
//        ocor_a *= ogamma;
//        ocor_b *= ogamma;
//        aluo = -ogg;

//        ref Face nface = ref N_faces[nfid];
//        double ngamma = N_Gamma.getFaceValue(nfid);
//        double ngg = ngamma * nface.gDiff;
//        double nfacecor = ngamma * nface.calNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_cells);
//        ncor_a *= ngamma;
//        ncor_b *= ngamma;
//        alun = -ngg;

//        ref double adiag = ref Local_A_Diag[iid];
//        ref double b = ref Local_b[iid];
//        adiag += ogg + ngg + ocor_a + ncor_a;
//        b += -(ocor_b + ncor_b);

//        mtool.r1.Set(-aluo, aluo, 0, ofacecor);
//        mtool.r2.Set(aluo, adiag, alun, b);
//        mtool.r3.Set(0, alun, -alun, nfacecor);
//        mtool.Eliminate();
//        O_A_LU_N[ofid] = -mtool.r1.a1;//用以计算A_Diag，还得取反
//        N_A_LU_N[nfid] = -mtool.r3.a3;
//        BetZone_A_LU_O[iid] = mtool.r1.a3;
//        BetZone_A_LU_N[iid] = mtool.r3.a1;
//        O_bFaceCommon[ofid] = mtool.r1.b;
//        N_bFaceCommon[nfid] = mtool.r3.b;

//        O_A_LU_O[ofid] = N_A_LU_O[nfid] = 1;//这是与对面区域的系数
//    }
//}
//O_Index.FaceIdx.GetSeg(out int ofid, out int endp1, taskIndex);
//int nfid = ofid + NFaceBaseO, iid = ofid + IFaceBaseO;
//for (; ofid < endp1; ofid++, nfid++, iid++)
//{
//    ref double alun = ref Local_A_LU_N[iid];
//    ref double aluo = ref Local_A_LU_O[iid];

//    ref Face oface = ref O_faces[ofid];
//    double ogamma = O_Gamma.getFaceValue(ofid);
//    double ogg = ogamma * oface.gDiff;
//    double ofacecor = ogamma * oface.calNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_cells);
//    ocor_a *= ogamma;
//    ocor_b *= ogamma;
//    aluo = -ogg;

//    ref Face nface = ref N_faces[nfid];
//    double ngamma = N_Gamma.getFaceValue(nfid);
//    double ngg = ngamma * nface.gDiff;
//    double nfacecor = ngamma * nface.calNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_cells);
//    ncor_a *= ngamma;
//    ncor_b *= ngamma;
//    alun = -ngg;

//    ref double adiag = ref O_A_Diag[oface.NeighborCellIndex];
//    ref double b = ref O_b[oface.NeighborCellIndex];
//    adiag += ogg + ngg + ocor_a + ncor_a;
//    b += -(ocor_b + ncor_b);

//    double commonTerm = 1 / adiag;
//    O_A_LU_O[ofid] = N_A_LU_O[nfid] = -aluo * alun * commonTerm;//这是与对面区域的系数
//    O_A_LU_N[ofid] = aluo + aluo * aluo * commonTerm;//用以计算A_Diag，还得取反
//    N_A_LU_N[nfid] = alun + alun * alun * commonTerm;
//    commonTerm *= b;
//    O_bFaceCommon[ofid] = ofacecor - aluo * commonTerm;
//    N_bFaceCommon[nfid] = nfacecor - alun * commonTerm;
//}
//TasksControl.AddDouble(ref N_A_Diag[nface.OwnerCellIndex], ngg);
//TasksControl.AddDouble(ref N_b[nface.OwnerCellIndex], nfacecor);

//O_A_Diag[oface.NeighborCellIndex] += ogg + ngg;
//O_b[oface.NeighborCellIndex] -= ofacecor + nfacecor;
//TasksControl.AddDouble(ref O_A_Diag[oface.OwnerCellIndex], ogg);
//TasksControl.AddDouble(ref O_b[oface.OwnerCellIndex], ofacecor);

//O_Index.FaceIdx.GetSeg(out int oface, out int endp1, taskIndex);
//for (; oface < endp1; oface++)
//{
//    int ocell = O_faces[oface].NeighborCellIndex;
//    int nface = oface + deltaFaceIndex;
//    int ncell = ocell + deltaCellIndex;

//    N_x[ncell] = O_x[ocell] = (O_b[ocell]
//        - O_A_LU[oface] * O_x[O_faces[oface].OwnerCellIndex]
//        - N_A_LU[nface] * N_x[N_faces[nface].OwnerCellIndex]
//          ) / O_A_Diag[ocell];
//}
//public void Ini_Ab(int taskIndex)
//{
//    if (S is Value1Base s1)
//    {
//        if (s1 is Value1Const Stemp && Stemp.number == 0)
//        {
//            O_Index.CellIdx.GetSeg(out int i1, out int endp1, taskIndex);
//            for (; i1 < endp1; i1++)
//            {
//                O_b[i1] = 0;
//                O_A_Diag[i1] = 0;
//                int nidx = i1 + deltaCellIndex;
//                N_b[nidx] = 0;
//                O_A_Diag[nidx] = 0;
//            }
//        }
//        else
//        {
//            O_Index.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//            for (; i1 < endp1; i1++)
//            {
//                ref Face face = ref O_faces[i1];
//                O_b[face.NeighborCellIndex] = s1.getFaceValue(i1) * face.Area;
//                O_A_Diag[face.NeighborCellIndex] = 0;
//                int nidx = face.NeighborCellIndex + deltaCellIndex;
//                N_b[nidx] = 0;
//                O_A_Diag[nidx] = 0;
//            }
//        }
//    }
//    else
//    {
//        throw new NotImplementedException();
//    }

//}

//public void Cal_Ab(int taskIndex)
//{
//    O_Index.FaceIdx.GetSeg(out int i1, out int endp1, taskIndex);
//    for (; i1 < endp1; i1++)
//    {
//        ref Face oface = ref O_faces[i1];
//        double ogamma = O_Gamma.getFaceValue(i1);
//        double ogg = ogamma * oface.gDiff;
//        // double ofacecor = ogamma * oface.calNonOrthogonalityCorrection_bound(O_grad, O_x);
//        double ofacecor = ogamma * oface.calNonOrthogonalityCorrection_bound(out double ocor_a, out double ocor_b, O_grad, O_x, O_cells);
//        ocor_a *= ogamma;
//        ocor_b *= ogamma;
//        O_A_LU[i1] = -ogg;
//        TasksControl.AddDouble(ref O_A_Diag[oface.OwnerCellIndex], ogg);
//        TasksControl.AddDouble(ref O_b[oface.OwnerCellIndex], ofacecor);

//        int nfidx = i1 + deltaFaceIndex;
//        ref Face nface = ref N_faces[nfidx];
//        double ngamma = N_Gamma.getFaceValue(nfidx);
//        double ngg = ngamma * nface.gDiff;
//        // double nfacecor = ngamma * nface.calNonOrthogonalityCorrection_bound(N_grad, N_x);
//        double nfacecor = ngamma * nface.calNonOrthogonalityCorrection_bound(out double ncor_a, out double ncor_b, N_grad, N_x, N_cells);
//        ncor_a *= ngamma;
//        ncor_b *= ngamma;
//        N_A_LU[nfidx] = -ngg;
//        TasksControl.AddDouble(ref N_A_Diag[nface.OwnerCellIndex], ngg);
//        TasksControl.AddDouble(ref N_b[nface.OwnerCellIndex], nfacecor);

//        //O_A_Diag[oface.NeighborCellIndex] += ogg + ngg;
//        //O_b[oface.NeighborCellIndex] -= ofacecor + nfacecor;
//        O_A_Diag[oface.NeighborCellIndex] += ogg + ngg + ocor_a + ncor_a;
//        O_b[oface.NeighborCellIndex] -= ocor_b + ncor_b;
//    }
//}

//protected void cal_x(int oface)
//{
//    int ocell = O_faces[oface].NeighborCellIndex;
//    int nface = oface + deltaFaceIndex;
//    int ncell = ocell + deltaCellIndex;

//    N_x[ncell] = O_x[ocell] = (O_b[ocell]
//        - O_A_LU[oface] * O_x[O_faces[oface].OwnerCellIndex]
//        - N_A_LU[nface] * N_x[N_faces[nface].OwnerCellIndex]
//          ) / O_A_Diag[ocell];
//}
//public void Set_A_Diag(int taskIndex)
//{
//    O_Index.getFaceSection(out int i1, out int endp1, taskIndex);
//    for (; i1 < endp1; i1++)
//    {
//        O_A_Diag[O_faces[i1].NeighborCellIndex] = -O_A_LU[i1] - N_A_LU[i1 + deltaFaceIndex];
//    }
//}
//public void Correct(int taskIndex)
//{
//    O_Index.getFaceSection(out int i1, out int endp1, taskIndex);
//    for (; i1 < endp1; i1++)
//    {
//        O_b[O_faces[i1].NeighborCellIndex] -= O_faceCommon[i1] + N_faceCommon[i1 + deltaFaceIndex];
//    }
//}
//public void Cal_x_Jacobi_Step1(int taskIndex)
//{
//    cal_x(x_owner, x_neighbor, x_OwnerCommon, x_NeighborCommon, taskIndex);
//}
//public void Cal_x_Jacobi_Step2(int taskIndex)
//{
//    cal_x(x_OwnerCommon, x_NeighborCommon, x_owner, x_neighbor, taskIndex);
//}


//protected void cal_x(double[] fromArray_O, double[] fromArray_N, double[] toArray_O, double[] toArray_N, int taskIndex)
//{
//    int oface, endp1;
//    OwnerIndex.getFaceSection(out oface, out endp1, taskIndex);
//    for (; oface < endp1; oface++)
//    {
//        int ocell = faces_owner[oface].NeighborCellIndex;
//        int nface = oface + deltaFaceIndex;
//        int ncell = ocell + deltaCellIndex;

//        toArray_N[ncell] = toArray_O[ocell] = (b_owner[ocell]
//            - A_LU_owner[oface] * fromArray_O[faces_owner[oface].OwnerCellIndex]
//            - A_LU_neighbor[nface] * fromArray_N[faces_neighbor[nface].OwnerCellIndex]
//              ) / A_Diag_owner[ocell];
//        //Console.WriteLine($"toArray{toArray_N[ncell]},aluo{A_LU_owner[oface]},alun{A_LU_neighbor[nface]},adiag{A_Diag_owner[ocell]},b{b_owner[ocell]},fromo{fromArray_O[faces_owner[oface].OwnerCellIndex]},formn{fromArray_N[faces_neighbor[nface].OwnerCellIndex]}");
//    }
//}
//this.O_x = ownerScalar.x;
//this.O_A_Diag = ownerScalar.A_Diag;
//this.O_A_LU = ownerScalar.A_LU_Diff;
//this.O_b = ownerScalar.b;
//this.O_Gamma = ownerScalar.Gamma;

//this.N_x = neighborScalar.x;
//this.N_A_Diag = neighborScalar.A_Diag;
//this.N_A_LU = neighborScalar.A_LU_Diff;
//this.N_b = neighborScalar.b;
//this.N_Gamma = neighborScalar.Gamma;
//set_b_owner_value(ref faces_owner[i1], i1);
//protected void set_b_owner_value(ref Face face, int faceIndex)
//{
//    b_owner[face.NeighborCellIndex] = S.getFaceValue(faceIndex) * face.Area;
//}
//public void pass_x_e_ToCoarser(out double[] e_ownerfine, out double[] e_neighborfine)
//{
//    e_ownerfine =this.x_owner;
//    e_neighborfine = this.x_neighbor;
//}
//- cells_owner[ocell].calALUProductValue(faces_owner, A_LU_owner, x_owner)
//- cells_neighbor[ocell + deltaIndexOfN_O].calALUProductValue(faces_neighbor, A_LU_neighbor, x_neighbor);
//public void pass_A_LU_owner_ToCoarser(out double A_lu_owner, ref CoarseFace ownerCoarser)
//{
//    double aluo = 0;
//    foreach (int fineface in ownerCoarser.FineFaces_O2N)
//    {
//        aluo += A_LU_owner[fineface];
//    }
//    A_lu_owner = aluo;
//}

//public void pass_A_LU_neighbor_ToCoarser(out double A_lu_neighbor, ref CoarseFace neighborCoarser)
//{
//    double alu = 0;
//    foreach (int fineface in neighborCoarser.FineFaces_O2N)
//    {
//        alu += A_LU_neighbor[fineface];
//    }
//    A_lu_neighbor = alu;
//}
//protected void set_b_value_owner(ValueBase value, IndexNumber indexNumber, int taskIndex)
//{
//    int i1, endp1;
//    indexNumber.getFaceInterval(out i1, out endp1, taskIndex);
//    for (; i1 < endp1; i1++)
//    {
//        b_owner[faces_owner[i1].NeighborCellIndex] = value.getFaceValue(i1) * faces_owner[i1].Normal.Length;
//    }
//}